home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
Almathera Ten Pack 3: CDPD 3
/
Almathera Ten on Ten - Disc 3: CDPD3.iso
/
fish
/
001-100
/
001-025
/
013
/
halley.bas
< prev
next >
Wrap
BASIC Source File
|
1995-03-17
|
6KB
|
213 lines
1 ' Comet ephemeris
2 gosub 10000
20 ' LET PI=3.1.159267
30 LET CO$="COMET HALLEY"
40 LET PH=1986.11
50 LET PL = 170.011
60 LET AN = 58.1453
70 LET PY = 76.0081
80 LET SM = 17.9435
90 LET EO = 0.967267
100 LET IO = 162.239
110 REM -------------------------
120 PRINT TAB (10);CO$
130 PRINT "- - - - - - - - - - - - - - -"
140 PRINT " EPHEMERIS FOR DATES"
150 PRINT " BETWEEN 1946 AND 2026"
160 PRINT "By Roger Browne"
165 Print "Adapted for Amiga by Tim Holloway":print
170 REM ========================
180 ' INPUT THE DATE
190 '
200 PRN=0:PRINT "Input year";
210 INPUT Y
220 IF Y<1946 OR Y > 2026 THEN 200
230 PRINT "Input month";
240 INPUT M
250 IF M < 1 OR M > 12 THEN 230
260 PRINT "Input day";
270 INPUT D
280 PRINT
290 PRINT INVERSE (1);" Calculating ... "
300 ' CALCULATIONS FOR THE COMET
310 '
320 LET X = PH
330 IF Y >= 1986 THEN Z=1984
340 IF Y < 1986 THEN Z = 1988
350 IF Y >= 1986 THEN S=0
360 IF Y < 1986 THEN S = 1
370 GOSUB 1780
380 DS = N
390 B = (360/PY)*(N/365.25)
400 K=B
410 GOSUB 1930
420 B=(K*PI)/180
430 E=B
440 Y1=EO
450 Q=E-(Y1*SIN(E))-B
460 IF ABS(Q) <= 1.7E-05 THEN 500
470 U=Q/(1-(Y1*COS(E)))
480 E=E-U
490 GOTO 450
500 V=(SQR((1+Y1)/(1-Y1))*(SIN(E/2)/COS(E/2)))
510 V=2*ATN(V)
520 V1=(V*180)/PI
530 L=V1+PL
540 R=SM*(1-(Y1*Y1))/(1+Y1*COS(V))
550 F=L-AN
560 F2=IO
570 F1=F*PI/180
580 F2=F2*PI/180
590 I=SIN(F1)*SIN(F2)
600 I=ATN(I/SQR(-I*I+1))
610 P=ATN((SIN(F1)/COS(F1))*COS(F2))
620 P1=P*180/PI+AN
630 IF F>= 90 AND F<= 270 THEN P1=P1+180
640 IF P1 < 0 THEN P1=P1+360
650 P=P1*PI/180
660 R2=R*COS(I)
670 '===========================
680 ' Calculations for the Earth
690 '
700 X=1975
710 IF Y>= X THEN Z=1972
720 IF Y < X THEN Z=1976
730 IF Y>= X THEN S=0
740 IF Y<1976 THEN S=1
750 GOSUB 1780
760 T=(360/365.25)*(N/1.00004)
770 K=T
780 GOSUB 1930
790 T=K
800 T1=T*PI/180
810 C=0.01672
820 J=T+(360/PI)*C*SIN(T1-0.051943)
830 J=J+99.5343
840 IF J>360 THEN J=J-360
850 IF J<0 THEN J=J+360
860 h=((j-102.51044)*pi)/180
870 r1 = (1-c*c)/(1+c*cos(h))
880 REM ------------------------
890 REM Compute Ecliptic coordinates
900 REM ------------------------
910 U1=((P1-J)*PI)/180
920 U2=((J-P1)*PI)/180
930 IF (R2<R1) THEN 990
940 Q1=(R1*SIN(U1))
950 Q1=Q1/(R2-(R1*COS(U1)))
960 Q1=ATN(Q1)
970 Q2=(Q1*180)/PI+P1
980 GOTO 1030
990 Q3=R2*SIN(U2)
1000 Q3=Q3/(R1-(R2*COS(U2)))
1010 Q3=ATN(Q3)
1020 Q2=(Q3*180)/PI+J+180
1030 IF Q2>360 THEN Q2=Q2-360
1040 IF Q2<0 THEN Q2=Q2+360
1050 Q4=Q2*PI/180
1060 Q5=(R2*(SIN(I)/COS(I))*SIN(Q4-P))
1070 Q5=Q5/(R1*SIN(U1))
1080 Q5=ATN(Q5)
1090 REM ------------------------
1100 REM Convert to equatorial coordinates
1110 REM ------------------------
1120 E1=0.40893064
1130 L1=SIN(Q5)*COS(E1)
1140 L1=L1+COS(Q5)*SIN(E1)*SIN(Q4)
1150 M1=ATN(L1/SQR(-L1*L1+1))
1160 Y2=M1*180/PI
1170 B1=(SIN(Q4)/COS(Q4))*COS(E1)
1180 B1=B1-(((SIN(Q5)/COS(Q5))*SIN(E1))/COS(Q4))
1190 G=ATN(B1)
1200 H1=G*180/PI
1210 I1=INT(Q2/90)
1220 J1=INT(H1/90)
1230 IF (I1-J1)=4 OR (I1-J1)=1 THEN H1=H1+360
1240 IF (I1-J1)=2 OR (I1-J1)=3 THEN H1=H1+180
1250 IF (I1-J1)=-4 THEN H1=H1+360
1260 IF (I1-J1)=-2 THEN H1=H1-180
1270 LET N1=H1/15
1280 W=INT((N1-INT(N1))*60+0.5)
1290 IF W=60 THEN N1=N1+1
1300 IF W=60 THEN W=0
1310 K1=ABS(Y2)
1320 W1=INT((K1-INT(K1))*60+0.5)
1325 G1=INT(K1)
1330 IF W1=60 THEN G1=G1+1
1340 IF W1=60 THEN W1=0
1360 IF (Y2<0) AND (G1<1) THEN W1=-W1
1370 D1=R1*R1+R2*R2
1380 D1=D1-(2*R1*R2*COS(U1))
1390 D2=SQR(D1)
1400 R3=D2/COS(I)
1410 K9=R
1420 GOSUB 2040
1430 R=K9
1440 K9=R3/10
1450 GOSUB 2040
1460 R3=K9*10
1470 MO=4.1:N=3.1
1480 IF DS<0 THEN MO=5:N=4.44
1490 MA=MO+5*0.4343*LOG(R3)
1500 MA=MA+N*2.5*0.4343*LOG(R)
1510 MA=(INT(10*MA))/10
1520 IF Y2<0 THEN G1=-G1
1530 REM ------------------------
1540 REM Print ephemeris for date
1550 REM ------------------------
1560 PRINT "--------------------"
1570 PRINT "Data for ";CO$:PRINT
1575 gshape (200,150),pic%
1580 print "Date: ";M;"/";D;"/";Y
1590 print "Days to perihelion: ";int(ds)
1600 print
1610 print "Coordinates: epoch 1950"
1620 print "RA ..... ";INT(N1);" hrs ";W; " min"
1630 print "Dec .... ";G1;" deg ";w1;" min":print
1650 print "Distances:"
1660 print "Comet to sun: . ";R;" A.U."
1670 print "Comet to Earth: ";R3;" A.U."
1680 print
1690 print "Predicted mag. ";MA
1700 print "--------------------"
1710 INPUT "Run another date";DATE$
1730 IF (DATE$<>"Y") AND (DATE$<> "y") THEN STOP
1740 GOTO 200
1770 REM ------------------------
1780 ' DAYS TO PERIHELION
1785 A=(Y-Z)/4
1790 A1=INT(A+S)
1800 N=365*(Y-X+S)+A1
1810 IF INT(A)<> A THEN 1830
1820 IF (M=2) AND (D<29) OR (M=1) THEN N=N-1
1830 IF M>2 THEN 1870
1840 M2=M-1
1850 M2=31*M2
1860 GOTO 1890
1870 M2=M+1
1880 M2=INT(30.6*M2)-63
1890 N=N+M2+D-365*S
1900 RETURN
1910 REM ------------------------
1920 REM PLACE BETWEEN 0 AND 360 DEG
1930 IF K<0 THEN 1950
1940 IF K>360 THEN 1980
1950 K=K+360
1960 IF K>=0 THEN 2010
1970 GOTO 1950
1980 K=K-360
1990 IF K<= 360 THEN 2010
2000 GOTO 1980
2010 RETURN
2020 REM ------------------------
2030 REM ROUND OFF SUBROUTINE
2040 K9=K9*1000
2050 K9=INT(K9+0.5)
2060 K9=K9/1000
2070 RETURN
10000 ' AMIGA GRAPHICS - OPTIONAL (TFH)
10005 scnclr
10010 dim pic%(138):bload "HALLEY.PIC",varptr(pic%(0))
10020 gshape (220,20),pic%
10030 return